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In this paper, we study shocks and related transitions in asymmetric simple exclusion processes 
of particles with nearest neighbor interactions. We consider two kinds of inter-particle interactions. 
In one case, the particle-hole symmetry is broken due to the interaction. In the other case, particles 
have an effective repulsion due to which the particle-current-density drops down near the half filling. 
These interacting particles move on an one dimensional lattice which is open at both the ends with 
injection of particles at one end and withdrawal of particles at the other. In addition to this, there are 
possibilities of attachments or detachments of particles to or from the lattice with certain rates. The 
hydrodynamic equation that involves the exact particle current-density of the particle conserving 
system and additional terms taking care of the attachment-detachment kinetics is studied using the 
techniques of boundary layer analysis. 



I. INTRODUCTION 

Asymmetric simple exclusion process (ASEP) com- 
prises of particles performing biased hopping in a pre- 
ferred direction on an one-dimensional chain. This bi- 
ased hopping gives rise to a finite particle current due to 
which the detailed balance is violated. In its simplest ver- 
sion, the particles respect mutual exclusion due to which 
a lattice site cannot be occupied by more than one par- 
ticle. Recently, ASEP with nonconservation of particles 
[l[ have drawn a lot of attention because these systems 
have better resemblance with the motors participating in 
biological transports inside the cell Q than the particle 
number conserving systems. The particles undergoing 
ASEP arc analogous to the molecular motors performing 
directed motion on tracks laid by the biopolymers [l], H[ • 
The conservation of particle number may be violated by 
absorption (or evaporation) of particles to (or from) the 
chain, a process, that is similar to the attachment or de- 
tachment of motors from the filament. 

One of the important issues involving these nonequi- 
librium processes on an open system is the boundary 
induced phase transition [4[. Boundary induced phase 
transitions are possible for bothparticle conserving Q 
and nonconscrving ASEP 0, i a a i with the latter 
one having a richer phase diagram with new phases. The 
other important aspect which influences the phase dia- 
gram is the interaction between the particles in addition 
to the mutual exclusion These interactions affect 

the particle current along the lattice. Since the parti- 
cle current carries the boundary information to the bulk, 
boundary induced phase transitions are likely to be af- 
fected by the change in the inter-particle interaction. 

The ASEP we consider here is defined on an one di- 
mensional lattice of size I with N lattice points [ll[ ■ Each 
site can be empty or occupied by a particle. A particle 
on site i can hop to site i + 1 with unit rate provided 
the target site is empty. A steady flow of particles along 
this one-dimensional channel is maintained by injection 
and withdrawal of particles at certain rates at left and 
right ends of the lattice respectively. The boundaries are 
coupled with particle reservoirs with fixed particle den- 



sities. We assume that the rates are such that the left 
and the right boundaries have particle densities a and 
7 respectively. The particle number in the bulk is not 
conserved since the particles are allowed to detach (or 
attach) from (or to) the chain with rates u>d K). In 
the steady state, depending on the values of a and 7, 
the system is found to exist in different phases that are 
characterized by the particle density p and the current- 
density j(p). These phases are represented in the phase- 
diagram in the space of a and 7. Monte Carlo simula- 
tions and mean field analysis of this model [l], Q show 
low-density, high-density and maximal current phases in 
addition to a phase, known as the shock phase, where the 
low-density and high-density regions coexist. The density 
profile across the lattice has a jump discontinuity (shock) 
from a low-density region to a high-density region. It has 
also been shown that the transition to the shock phase 
can become critical under special values of the boundary 
parameters [l[ . Shocks (domain walls) separating differ- 
ent steady-state density profiles have been observed in 
ASEP without the absorption-desorption kinetics. How- 
ever, in this case, the position of a shock on the lattice 
fluctuates and they can, in general, move on the lattice. 
The dynamics of such shocks [HI E2| drive the system to 
settle in one of its steady states which are known as high- 
density, low-density and maximal current phases (in case 
of particles having only mutual exclusion) . Studying the 
dynamics of such shocks helps in understanding the na- 
ture of the phase transitions in such systems. Unlike this 
particle conserving ASEP, here the shocks are localized 
@ and there are stationary states in the phase dia- 
gram where the density profiles have such shocks. The 
shock is upward in the sense that the low-density phase 
appears on the left of the shock and the high-density 
phase appears on the right of the shock. The position 
of these stationary shocks is determined by a, 7 and the 
absorption-desorption rates. A boundary layer analysis 
shows that the transition to the shock phase from a low- 
density phase can happen through a critical deconfrne- 
ment of a boundary layer from the edge of the system 
0- In addition, it has been shown in a general way that 
this transition to the shock phase is associated with a 
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dual transition [8|, |9( . The dual transition is a boundary 
transition across which the slope of the boundary layer 
changes sign with the shape of the bulk density profile 
remaining unchanged. These boundary and shock transi- 
tion lines meet at the critical point. The presence of the 
dual transition helps in developing a general approach 
through which one can classify the phase-diagrams into 
different categories. This has been done in reference @ 
using the knowledge of zeros of a set of coarse grained 
functions present in the equation describing the steady 
state density profile. 

The effect of the inter-particle interactions on the par- 
ticle conserving ASEP has been found to be quite drastic 
[n|. Typically one considers interactions of the follow- 
ing nature [1 31 ] . The particles hop to the right with bulk 
hopping rates 

10 0^0010 with rate 1 + 6 (1) 

1100^1010 with rate 1 + A (2) 

0101^0011 with rate 1 - A (3) 

1101-^1011 with rate 1-8. (4) 

Here, 1 and imply an occupied and an unoccupied site 
respectively A = 6 = 0, is the case of ASEP with 
only mutual exclusion. In this case, the current-density 
j(fi) = p(l — P) respects particlc-holc symmetry and has 
a maximum at p = 1/2. An exact solution of the particle 
conserving problem [|| reveals that there are primarily 
three phases known as the low-density, high-density and 
the maximal current phases. The case of A — > 1 im- 
plies a repulsion among the particles due to which the 
current versus density plot is expected to show a min- 
imum near the half-filling. In the range < A < 1, 
as A increases, the current-density plot develops double 
maxima from a single maximum picture. With a mini- 
mum in the current versus density plot, the phase dia- 
gram contains a minimum-current phase [l(| in addition 
to the low-density, high-density and the maximal current 
phases. At this minimum current phase, the bulk den- 
sity has a value that corresponds to the minimum value 
of the current. The maximal-current phase, in general, 
is also expected to be of two different types in which 
the bulk density corresponds to two different maxima of 
the current. Equations [1] and 2] involving 8 imply that 
the particle current is the same as the vacancy current 
for 8 = 0. A finite 8 breaks this equality since 8 > 
and 8 < imply larger particle or vacancy currents re- 
spectively. This particle-hole asymmetry introduces an 
asymmetry in the current versus density plot. 

Since the particle nonconscrvation and the inter- 
particle interaction both individually leave significant im- 
pact on ASEP, their combined effect is expected to be 
interesting. Although the phase diagram of the parti- 
cle nonconserving version of the interacting ASEP has 
not been probed, numerical simulations exhibit some new 
features such as downward shock (see figure [T|) and dou- 
ble shocks [HI in the density profile. In case of down- 
ward shock, the shock is a discontinuity in the density 
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FIG. 1: A diagram for the density profile with a downward 
shcok. 



profile separating a high-density part on the left and a 
low-density part on the right. The density profile with a 
double shock consists of two successive upward jumps 
separating a low-density part on the left and a high- 
density part on the right. The particular feature that 
interests us here is the downward shock in the density 
profile. 

In this paper, we look into the details of the shocks for 
two limiting cases (i) A = and 8 = 1 and (ii) A = 1 and 
6 = using a boundary layer analysis. These problems 
serve as the limiting cases of the most general interact- 
ing case which is technically more difficult to handle. The 
first case is useful for the verification of our earlier pre- 
dictions related to shockening transition, dual boundary 
transition and criticality. We show that although the 
shape of the density profile is different in this case, dif- 
ferent phase transitions and the phase diagram are quali- 
tatively same as those of only mutual exclusion case with 
Lj a ^ LOd- Further, the asymmetry in the current versus 
density plot does not alter the exponents which remain 
the same as those of the pure mutual exclusion case with 
Lj a 7^ cod- In the second case, our focus is particularly 
on the shape of the density profile that has a downward 
shock, location of the shock and specific conditions on 
the boundary parameters a and 7 for which such shock 
can be observed. The aim is also to show how the earlier 
analysis for the upward shock can be appropriately modi- 
fied to understand the details of the downward shock. We 
show that the earlier predictions [l4| about the ranges in 
which the values of the densities at the upper and lower 
end of the shock must lie, are satisfied naturally in the 
boundary layer analysis. Although the boundary layer 
analysis presented here can be used to obtain other fea- 
tures such as phases with upward shocks, entire phase di- 
agram with the details of the phase boundaries, we plan 
to present those with supportive numerical analysis of 
the system in a later publication. 

The paper is organized as follows. In the next section, 
we discuss some of the details of the model. In section 
III, we consider the first case, namely that of A = 
and 8 = 1. Some results related to the density profiles 
and phase diagrams for this case are presented in the 
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Appendix. In section IV, we consider the case of A = 1 
and 6 = 0. We conclude the paper with a discussion in 
section V. 



II. MODEL 

In general, for a particle conserving system, one can 
write down a discrete continuity equation 



dpk 
dt 



Jk-i - Jk, 



(5) 



where jk is the average current across the bond between 
the fcth and k — 1th sites and pk is the average parti- 
cle occupancy at the fcth site. The average current takes 
care of the particle hopping rule that is problem specific. 
In the presence of particle absorption-desorption kinet- 
ics, we may supplement this equation with other terms 
responsible for particle absorption and desorption. One 
can obtain a simple, solvable continuum model from the 
resulting equation by making a mean field approximation 
that neglects particle correlations. In case of only mutual 
exclusion, the current-density obtained from the mean- 
field approach is the exact one and as a consequence of 
this there is a good agreement in results obtained from 
the mean-field analysis and the Monte-Carlo simulation. 
The mean-field approach, however, does not produce cor- 
rect results in the presence of additional interactions con- 
sidered here. This is because the mean-field approach 
fails to produce even the basic qualitative features of the 
current-density. 

In view of these issues, we follow the approach of hy- 
drodynamic equation [14j . In this approach, we start 
with the thermodynamic limit, N — > oo with the lat- 
tice spacing a = l/N — > 0. For the entire analysis in 
the following, we choose 1 = 1. In the continuum limit, 
k — > x = ka and lattice — > t = t latticc a, the average 
particle density p(x, t) satisfies the continuous version of 
the continuity equation. With the particle nonconserving 
terms as additional source and sink terms, we have the 
continuous equation 

Here il a = uo a N and ild = ^d,N. We use the exact 
stationary current in the hydrodynamic equation with 
the assumption that the bulk has sufficient time to relax 
between the particle absorption and desorption events. 
This is possible when the absorption and desorption hap- 
pen at a very low rate. 

The current-density for the interacting system has 
been derived in the past [lj| by using the fact that the 
stationary distribution is given by the equilibrium distri- 
bution of an one dimensional Ising model. The average 
current across the bond between k and k + 1th sites is 
given by 

jk = (1 + <5)<0100) + (1 + A)(1100) + 
(1-A)(0101) + (1-<J)(1101), (7) 



where the four sites in the equation are the k—1, k, k + 
1, k + 2th sites. The four-point averaging needs to be 
done with respect to the stationary measure. Since the 
details of the calculations are available in [l(| , we quote 
the final result for the current-density. In the thermody- 
namic limit, N — > oo, the current-density is 



. A[l + ^(l-2p)]-A A /4p(l-p) 
3 = . ( g ) 



where 



+ (: 



y/4p(l ~ p) Ml-P) 



- 1 



a/2 



1 



(9) 



For the mutual exclusion case, A = 5 = 0, one obtains 
the exact current-density j = p(l — p). This is also the 
current-density one finds through a mean-field approach 
starting with jk = (10). Mean-field amounts to ignoring 
the correlation and assuming jk = (1)(0). For 6=1 
and A = 0, the current-density is j = 2p(l — p) 2 . The 
current-density vanishes in the completely unoccupied, 
p — > 0, and completely occupied, p — ► 1, limits with 
a maximum at p = 1/3. The lack of the particle-hole 
symmetry is reflected in the asymmetric nature of the 
current-density around p = 1/2. For A = 1 and 6 = 0, 
the current-density is 



j = v(i-v)/{i + v), 



where 



y 



y/l - 4p(l - p). 



(10) 



(11) 



The current-density is symmetric about p = 1/2, and has 
two maxima at p = .707 and p = .293 (see figure [2|). 




FIG. 2: The plot of the particle current versus particle density 
for 6 = and A = 1. 

Equation |6]) is subjected to the boundary conditions 
p(x = 0) = a and p(x = 1) = 7. Although the first order 
equation obtained from equation © under the steady- 
state condition, dp(x,t)/dt = 0, describes the shape of 
the steady-state density profile as a function of x, this 
equation, in general, does not have a smooth solution 
satisfying two boundary conditions. This problem can 
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be overcome by adding a second order term to this equa- 
tion. This second order term is introduced with a pref- 
actor which becomes vanishingly small in the thermody- 
namic limit. It can be shown that such an infinitesimal 
second order term is generated anyway as one takes the 
continuum limit of (|5|) and retains terms up to O(a). The 
second order term ensures a smooth solution for the hy- 
drodynamic equation through the formation of shocks or 
boundary layers over a region of width of the order of 
the lattice spacing, 1/N. In the steady-state, the final 
form of the equation describing the shape of the density 
profile is, therefore, 



e -^-^ + n a {i- P {x))-n dP {x) = Q, (12) 



where e is a factor of the 0(1 /N) and e — > in the 
thermodynamic limit, N — > oo. In the presence of the 
second order term, one can define an effective current as 



e dp 
2dx 



+ 3- 



(13) 



We discuss later the important role J plays in predicting 
the nature of the density profile. 

A knowledge of certain special densities obtained from 
equation (12) can be useful for us. In the presence of 
only absorption-desorption kinetics, with K = Q a /fld, 
the system settles into a stationary density, p = pl = 
K/(K +1) obtained from the zero of the particle non- 
conserving terms. It is possible for the system to settle 
into this constant equilibrium density profile in the pres- 
ence of hopping dynamics also provided the equation of 
motion allows such a density profile. This is what hap- 
pens in case of ASEP with only mutual exclusion. In 
this case, the current-density j(p) = p(l — p) has a max- 
imum (j'(p) \ Pm = 0)at p = 1/2 and with K = 1, one has 
Pl = Pm- This is a very special situation which gives rise 
to a phase in which the density profile can maintain the 
maximum current in the system. Thus for K — 1 , in ad- 
dition to the low-density, high-density and shock phases 
one also has a maximal current phase and coexistence of 
the maximal current phase with other phases. Different 
phases and different exponents associated with the phase 
transitions in the purely mutual exclusion case indicate 
the existence of a different universality class for K = 1. 
A recent paper Q shows that instead of the value of K, 
it is the relative magnitude of ph and p m that determines 
the phase diagram of the system. While pl = Pm ( for 
purely mutual exclusion case, this condition is fulfilled 
if K = 1) forms a special case where the maximal cur- 
rent phase appears, the phase diagrams for pl < p m and 
Pl > Pm are related through the particle-hole symmetry. 
In our entire analysis, we consider K = 1 (fi = CI4 = O). 
This, however, does not restrict us to any special univer- 
sality class as we maintain pl ^ p m . 



III. BASIC PRINCIPLES OF BOUNDARY 
LAYER ANALYSIS 



We mention here a few basic principles that have been 
employed in the following boundary layer analysis. Wc 
start from the basic hydrodynamic equation (|12p . Since 
the contribution of the second order term in this equation 
is very small in the thermodynamic limit, it is expected 
that the major part of the density profile is described 
by the solution of the first order equation obtained from 
equation (|12[) by ignoring the second order term. This 
solution, known as the outer solution, cannot satisfy both 
the boundary conditions in general. In order to satisfy 
the boundary conditions appropriately, there appear spe- 
cial narrow regions (of width of 0(e)) such as boundary 
layers or shock. In order to see these regions, it is nec- 
essary to rescale the position variable of the differential 
equation appropriately. These regions are described by 
inner solutions which are solutions of the rescaled equa- 
tion. Upon rcscaling equation Q12[) . the nonconservative 
terms acquire a prefactor of 0(e), due to which these 
terms are neglected. The inner solution is, therefore, de- 
termined by the first two terms of equation (fT2"|) . Differ- 
ent unknown constants, present in the solutions in differ- 
ent regions, are determined from the boundary conditions 
or by the conditions required for the smooth joining of the 
outer and inner solutions. This is the general scheme 15 1 
through which one can obtain an uniform approximation 
of the solution of equation (|12[) in the thermodynamic 
limit. 

Based on these principles, a few predictions about the 
presence or absence of the shock in the density profile 
can be made. The knowledge about the shape of the 
density profile with shock or the location of the shock, 
however, requires explicit calculations. Since the par- 
ticle nonconserving terms of equation (|12[) are not im- 
portant for the boundary layer or the shock region, wc 
expect the total current J to remain constant in this re- 
gion. For a shock separating two densities p\ and p ro , 
we have the total current J — j(p\ ) and J = j(p ro ) at 
the two edges of the shock and the constancy of the to- 
tal current J across the shock, demands j{p\ ) — j(Pm)- 
This condition is implemented explicitly in the boundary 
layer analysis in finding out the inner solutions describ- 
ing the shocks or boundary layers. As a consequence of 
this constancy condition, the path representing the shock 
on the current-density plot is horizontal. For an upward 
shock in the density profile (dp/dx > 0), one requires 
J < j(p(x)). The horizontal path representing the up- 
ward shock on the current-density plot has to lie below 
the curve j(p(x)). Therefore, we require at least one max- 
imum in the current and density relation to fulfill these 
conditions for an upward shock. Similarly, we expect 
j(p(x)) < J in case of a downward shock (dp/dx < 0). 
This requires a concave region in the current versus den- 
sity plot. In our case, due to the particle-hole symmetry, 
the center of the downward shock lies at p = 1/2. As a 
result, p\o and p TO are bounded as .5 < p\ Q < .707 and 
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.293 < p 10 < .5. In the boundary layer analysis, these 
conditions appear naturally as requirements for the sat- 
uration of the shock solution (inner solution) to the bulk 
on both sides of the shock. 



IV. BOUNDARY LAYER ANALYSIS FOR 8 = 1 

AND A = 

We start our analysis by finding out the outer solu- 
tions. The outer solution is the solution of equation 

-2(1-4/9 + 3^)^ + ^(1-2^=0, (14) 
dx 



obtained from equation (|12[) by ignoring the second order 
term. The solution is given in terms of the transcendental 
equation 



ff(Pout) = + c, (15) 



where 



9(p) = - g ( V - 5P) + 4 lo g[ 2 P - 1] ( 16 ) 

and c is a constant that can be found out from the bound- 
ary condition that the outer solution satisfies. As in the 
pure mutual exclusion case, there can be a phase where 
the outer solution satisfies the left boundary condition 
p(x = 0) = a. We consider this situation in the follow- 
ing. In this case, the outer solution is given by 



g(p) = nx + g(a). 



(17) 



The inner solution describing the boundary layers or 
shocks can be found from equation (|12p . by expressing it 
in terms of x = (x — xo)/e, where Xo represents the lo- 
cation of the inner solution. The nonconservative terms 
in the rescaled equation are negligible in the e — > limit. 
The inner solution is, therefore, the solution of the equa- 
tion 



1 dp 

2di 



2(p-2p 2 



P 3 ) + c, 



(18) 



where C is a constant. The saturation of the inner so- 
lution to p a = p ut{x — > 1) as x — > — oois ensured by 
choosing 



C=-2{ Po -2pl + pl). 
Equation (fT8|) can be rewritten as 
dp 



dx 



~ Po)(P ~ Pl)(p ~ P2), 



(19) 



(20) 



where 



pi,2 = ^[2- Po ±{ipo-^p 2 } 1/2 }. (21) 



The general solution of equation (|20"|) is 



log[ 



(p - Po) ipl -" 2 Hp ~ Pl) ( " a - po) 1 



(P-P 2 )(P1-Po) 

4(Po - pi)(po - P2)(pi - p 2 ){x + £), (22) 

where £ is a constant. There are possibilities of saturation 
of the inner solution to densities pi(p ) or P2(p ) as a; — > 
oo. These two densities are functions of a and Q through 
p . The possibility of saturation to p\ can be ruled out 
since for all values of p a in the range {0, 1}, pi > 1. The 
inner solution can, however, saturate to pi which remains 
in the realistic range (< 1) in the entire range of p a . In 
case of saturation, the approach to pi is given by 

P - P2+ (P2 - Po) {Pl - p2)/{Pl - p ° ] {p2 - Pl )^-P^/^-P^ 



exp[-4(p 2 - p )(pi - p 2 )(x + 0]- 



(23) 



Two length scales appear in the inner solution. The 
length scale £ describes the center of the inner solution. 
Using the boundary condition, p(x = 0) = 7, that equa- 
tion ([22]) must satisfy, we have 



1 



4(Po - Pl)(Pa - P2)(pl - P2) 

h-Po) pl - p2 h- P i) p2 - p ° 



log 



(24) 



The other length scale, which we denote as w, describes 
the approach of the inner solution to the saturation value. 
This length scale is 



w = VWPi - P2X/J2 - Po)]. 



(25) 



This case of outer and inner solutions satisfying the 
boundary condition at x = and x = 1 respectively 
is possible for low values of a. This phase, with a domi- 
nated bulk density for low values of a is the low-density 
phase. 

The boundary layer at x = 1 is unable to satisfy the 
right boundary condition if 7 > P2 and as a consequence 
of this, the boundary layer deconfines from the bound- 
ary at x = 1 as 7 > p2- To satisfy the right bound- 
ary condition, an outer solution joined smoothly to the 
boundary layer at the left and satisfying the right bound- 
ary condition at its right appears. The deconfincment of 
the boundary layer from the boundary as 7 exceeds p 2 
is the shockening transition that has been seen earlier 
in the pure mutual exclusion case. The phase bound- 
ary between the low-density and the shock phase can be 
determined from the condition 



p 2 («) = 7- 



(26) 



As the shock phase boundary is approached the length 
scale £ diverges logarithmically as 



f ~ log(7 - P2) 



(27) 



Equation (|20[) also exhibits a boundary transition that 
has been first identified as a dual transition in the pure 
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mutual exclusion case [8j and is expected to be present 
in general whenever there is a bulk shockening transi- 
tion through a deconfinemcnt of a boundary layer. This 
boundary transitions happens in the low-density phase 
and across this boundary transition, the slope of the 
boundary layer (inner solution) changes sign. It can be 
seen from equation (|20p , that the boundary layer has a 
positive slope if p a < 7 < pi and has a negative slope if 
7 < p Q . Thus, this change in the slope of the boundary 
layer happens across the boundary transition line 

Po{a) = 7. (28) 

We call the low-density phase with boundary layer hav- 
ing positive slope as low-density(l). The other part of the 
low-density phase is called as low-density (2). The bound- 
ary transition supports the duality theorem of reference 
Q that for every a if there is a bulk phase transition at 
7 = pii there is a boundary transition at 7 = p {a). As 
the boundary transition line is approached from cither 
side, the length scale £ diverges logarithmically as 



£ ~ log(po - 7)- 



The bulk shockening transition line and the boundary 
transition line meet at a point where 

Po = P2 = 7- (30) 

This is the critical point at which the length scale w 
diverges as 



1 



(P2 ~ Po) ' 



(31) 



Equation (|30)) . along with equation (|2T|) leads to the crit- 
ical value for 7 C = 1/3. The corresponding critical value 
of a can be determined by implementing the condition 
Po(c*e) = 1/3 in equation (JT7J). This leads to the following 
equation 



g(l/3) = Sl+g(a c ), 



(32) 



from which the O dependent a c can be determined. 

The analysis presented here agrees with the predic- 
tions of our earlier work. Although numerical values of 
a and 7 for the critical point or for other transitions are 
now shifted, the qualitative features of the phase dia- 
grams and phase boundaries remain same as that of the 
pure mutual exclusion case with Q a ^ 0^. We have 
shown here that the nature of divergences of different 
length scales associated with the shock transition and 
the boundary transition are same as those in the case of 
ASEP with only mutual exclusion 0, [1] . Since the cal- 
culation of other exponents are similar to those of ASEP 
with only mutual exclusion, we refer the reader to refer- 
ences @, HJ f° r this purpose. 



V. BOUNDARY LAYER ANALYSIS FOR 5 = 

AND A = 1 

For 5 = and A = 1, we look for those solutions of 
equation (|12p which support a downward shock in the 



density profile. In this system, the downward shock is 
seen due to the presence of a concave region in the current 
versus density plot. Due to the particle-hole symmetry, 
we expect the downward shock to be centered around 
p = 1/2. Further, there is a discontinuity in dj /dp at p = 
1/2, because of which we need to distinguish two regions, 
p > 1/2 and p < 1/2. As a consequence, the boundary 
layer analysis for the density profile differs significantly 
from that of the asymmetric case considered in the first 
:>art of this paper and also from our previous studies in 



part c 
0,1 



The bulk solutions on the left and right of the down- 
ward shock are described by the left and right outer so- 
lutions obtained below. The outer solutions can be ob- 
tained by solving equation (|12p with its first term ig- 
nored. The left outer solution, for which p > 1/2, is the 
solution of 



2(1 - 2y - y 2 ) dp 



(1 + y) 2 



dx 



0(1 - 2p) = 0, 



(29) with y = 2p — 1. The solution of the equation is 



9i(p 



lout 



= fix + C\ , 



where 



9i(p) = — + 21ogp-log(2p-l). 
P 



(33) 



(34) 



(35) 



ci can be fixed from the boundary condition, p(x = 0) = 
a, since the left outer solution satisfies the left boundary 
condition. 

The equation for the right outer solution can be found 
out from (|12[) in a similar way with the substitution 
y = 1 — 2/9. The right outer solution on the right of 
the downward shock is 



-(prout) = Ox + c 2 , 



(36) 



where 



9r(p) = 7 - 21og(p - 1) + log(2p - 1). (37) 

P- 1 

As before, C2 can be determined from the condition p(x = 
1) = 7. The plot of these two solutions for given a and 
7 is shown below. 

In the following, we find out the inner solutions that 
describe the downward shock. From the symmetry of the 
current-density plot about p = 1/2, we expect the center 
of the inner solution to correspond to p = 1/2. The in- 
ner solution that approaches the left outer solution can 
be obtained from equation (fl"2|) with y = 2p — 1, after 
expressing it in terms of x = (x — x s )/e, where x s repre- 
sents the center of the inner solution. With the rescaling 
of x, the O dependent term in equation (|12p drops out 
for having negligible contribution. The equation that de- 
termines the left inner solution is 



dp 
dx 



2 + 4p 2 



■di. 



(38) 



7 



P 



. 6 
. 55 

. 5 
.45 




FIG. 3: The plot of the left (upper curve) and right (lower 
curve) outer solutions. The left outer solution satisfies the left 
boundary condition a = .64. The right outer solution satisfies 
the right boundary condition 7 = .36. We have plotted the 
solutions till the other end of the lattice. 

We expect the inner solution to saturate to the left outer 
solution pi = p\ m t(x ~~ > x s —) as x — » —00. Therefore, 
we choose 

dl = 2 -±^. (39) 

Plo 

The solution of (|38|) is given by the transcendental equa- 
tion 

- y log (P - «°) + \ lo s( 2 » - 1) = (5 + c)(2pg, - 1J40) 

The constant c can be found out by demanding p(x = 
0) = 1/2. The inner solution should approach p\ as 
x — > — 00. This is possible if 2pf Q — 1 < 0. Thus on the 
high-density side, the shock should saturate to a density 
pio which is bounded as 

.707 > p lo > .5. (41) 

The right inner solution that approaches the right outer 
solution p ro = p rou t(x — > x s +) as x — + 00 can be found 
out in a similar way. Substituting y — 1 — 2p in equation 
(|12[) . we find the right inner solution 

^ log(p - Pro) - \ log(2pp IO - 1) = (2j& - + d)&42) 

where p = 1 — p. The constant d can be found out by 
demanding that p = 1/2 at x = 0. This solution also 
saturates to p ro as i — > 00, with the condition that 2p;? — 
1 < 0. Thus on the low-density side, the shock saturates 
to a value bounded as 

.29 < p ro < .5. (43) 

Conditions in (flTj) and (|13"|) are same as those predicted 
for a downward shock for the particle conserving inter- 
acting systems. The boundary layer analysis presented 
here shows that these conditions are necessary for the 
saturation of the inner solutions. 



In the fully symmetric situation (<5 = 0), the inner 
solution is centered at p = 1/2. This imposes another 
constraint 

Plo = 1 - Pro- (44) 

This condition also guarantees the matching of the slopes 
of the left and right inner solutions at x = x s . Equation 
([44)) allows us to find out the shock positions and the 
dependence of pi and p ro on a and 7. If the shock is 
formed at x — x s , we have 

5i (plo) = Clx s + Pi (a) (45) 
ffr(Pro) = Clx s + g T (j) - CI. (46) 

Using equations P3"|) . P6"| and (JUJ), we have the final 
equation, 

ffi(pio) = <?r(l - Pio) + n + \gi(a) - flr(7)], (47) 

that determines p\ a . Equation (|47|) and (|4T|) together set 
the condition for the formation of the shock. Knowing 
the value of pi from equation (|47"|) for given values of 
a, (3 and O, one can find out the height of the shock 

-Hshock = Plo — Pro = 2pi — 1. 

VI. SUMMARY 

In this paper, we present a boundary layer analysis for 
studying the shocks and associated transitions for asym- 
metric simple exclusion processes of interacting particles. 
The particles that hop in a particular direction on a fi- 
nite one-dimensional lattice interact mutually. These in- 
teractions are present in addition to the usual mutual 
exclusion among the particles. We consider two different 
kinds of interactions. In one case, particles have mutual 
repulsion. In the other case, particle-hole symmetry is 
broken by the interaction. We consider these two cases 
separately since in that case simple analytical solutions 
for the density profile with shock can be obtained. In ad- 
dition to these interactions, there arc possibilities of at- 
tachmcnt(detachmcnt) of particles to(from) the lattice. 
A steady flow of particles on the lattice is maintained 
through the injection of particles at one end and with- 
drawal of particles at the other end at certain rates. We 
assume that the rates are such that the particle densities 
at the two ends remain a and 7. Depending on the values 
of these boundary densities a and 7, the system, in the 
steady state, can be in different phases. These phases are 
characterized by the shapes of the particle density profiles 
across the lattice and the particle-current-densities. The 
hydrodynamic equation that involves the exact current- 
density of the particle conserving system is supplemented 
with additional particle absorption and desorption terms. 
In our entire analysis, we choose equal particle absorption 
and desorption rates. The resulting equation is, then, 
studied using the techniques of boundary layer analysis 
to obtain the steady-state density profiles. 
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In the case where the interaction breaks the particle- 
hole symmetry, the phase diagram and the phase tran- 
sitions arc qualitatively similar to that of the purely 
mutual exclusion case with unequal particle absorption- 
desorption rates. In our interacting system, the particle 
density profile shows a jump discontinuity or shock from 
a low value to a high value over an extended region on 
the a-7 space. As in the case of ASEP with only mu- 
tual exclusion, the transition to the shock phase happens 
through the deconfinement of the boundary layer. Asso- 
ciated to this transition to the shock phase, there is also 
a dual boundary transition across which the slope of the 
boundary layer changes sign with the bulk density pro- 
file remaining the same. The exponents characterizing 
various phase transitions are same as those of the ASEP 
with only mutual exclusion of particles. 

In case of repulsion among the particles, it is known 
that there exists a downward shock for certain values of 
the boundary densities. We obtain the analytical form of 
the density profile that has a downward shock. The so- 
lutions describing the shock region in the density profile 
saturate to the bulk part of the density profile exponen- 
tially. This exponential approach is possible if the bulk 
density at the edges of the shock remain within certain 
range. We further obtain conditions on the boundary pa- 
rameters for having such a downward shock, and also the 
location and the height of the shock. This work opens 
up the scope of analyzing more general interacting prob- 
lems. Even in the absence of any compact closed form 
analytical solutions for the density profile, the basic prin- 
ciples used in this problem should be applicable in more 
general cases. 



APPENDIX A: DENSITY PROFILES AT 
VARIOUS PHASES AND PHASE DIAGRAM FOR 

5 = 1, A = 

In this appendix, we present a few representative plots 
of the density profile at different phases. The density 
profile in the low-density phase (low-density (1)) has the 
shape as given in figured! Figure [5] shows the deconfincd 
boundary layer as 7 is increased from the shock transition 
value 7 = .43.. for a = .15. 

The deconfinement of the boundary layer to form the 
shock happens only for 7 > j c . For 7 < 7 C , as p a becomes 
equal to 1/3 with the increase in a, the shock starts form- 
ing with an outer solution appearing on the right of the 
shock and satisfying an effective right boundary condi- 
tion p(x = 1) — 7 C . The true boundary condition at 
x = 1 is satisfied by a decaying boundary layer simi- 
lar to the one present in the low-density phase. Since 
for all 7 < 7c, the effective boundary condition on the 
right for shock formation remains same, the low-density- 
shock phase boundary is vertical at a = a c for 7 < 1/3. 
The shock height increases from zero continuously as one 
enters the shock phase. The density profile typically ap- 
pears as in figure [HI 



0.5 r 

0.45 - 



0.4 
0.35 
Q. 0.3 




0.1 1 1 1 1 1 

0.2 0.4 0.6 0.8 1 

X 

FIG. 4: Plot of the density profile for a = .15, 7 = .42, 
Q, = .1, e — .03. At these values of the parameters, the system 
is in the low-density phase, just below shockening transition 
line. 



Q- 0.3 




0.2 0.4 0.6 0.8 1 

X 



FIG. 5: Plot of the density profile for a = .15, 7 = .45, 
Q = .1, and e = .03. At these values of the parameters, the 
boundary layer is just deconfined from the x — 1 boundary. 

The shock, after being formed at x — 1 , moves toward 
the bulk of the system as a is increased for a given 7 > j c . 
The shock continues to exist till the other shock phase 
boundary is reached. At the shock phase boundary, the 
outer solution satisfying the boundary condition p{x = 
1) = 7 spans almost the entire lattice except for a narrow 
region for the inner solution whose saturation value at 
x = just matches the boundary density a. With further 
increase of a, the shock disappears from the system and 
the system enters into a different phase where there is no 
shock. The density profile in this phase typically appears 
as in figure [7l Following the case of ASEP with only 
mutual exclusion, we call this phase as the high-density 
phase. In other words, as we approach the shock phase 
from the high-density side by decreasing a, the boundary 
layer deconfines from the x = boundary in the form of 
a shock as soon as a becomes smaller than the saturation 
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0.38 




X 



FIG. 6: Plot of the density profile for a = .20, 7 = .30, 
ft = .1, e = .009. At these values of parameters, there is a 
shock in the density profile but it is not formed through the 
deconfinement of the boundary layer. Two branches of outer 
solutions are joined discontinuously through the shock whose 
height increases as a increases with fixed 7. 




' 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 



X 

FIG. 7: Plot of the density profile for a = .1, 7 = .8, fl = .1 
and e = .034. At these values of parameters, the system is in 
the high-density phase. 

value of the inner solution at x = 0. As a result the phase 
boundary between the shock phase and the high-density 



phase is given by 

M) = a, (Al) 

where p' Q is the value of the outer solution at the bound- 
ary x = 0. Since the outer solution now satisfies the 
boundary condition at x = 1, p' a is a function of 7 and fl. 
Same equation is true for the phase boundary for 7 < j c 
except for the fact that here p' Q is independent of 7 since 
p(x = 1) = 7c. This leads to a vertical phase boundary 
between the shock and the high-density phases. With the 
phase boundaries specified as above, the phase diagram 
appears as in figure [5J The boundary transition exists 
in the high-density phase also. Depending on whether 
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FIG. 8: Quantitative plot of the phase diagram for 5=1 
and A = and fl — .1. Low-density (1) and low-density (2) 
are the two low-density phases in which the bulk profiles re- 
main same but the surface layers have different slopes. These 
two phases are thus separated by the boundary transition 
line. The phase boundary between the low-density and shock 
phases and the boundary transition line meet at the critical 
point {a = .183.., 7 = 1/3}. 

the value of a is larger or smaller than p' , the slope of 
the boundary layer at x = changes sign. Since the 
boundary transition lines in the high-density phase can 
be calculated following the principles mentioned above 
and also from reference Q, we skip those calculations 
here. 
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